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Abstract. 

The  purpose  of  this  research  is  to  study  the  effects  of 
material  properties,  thickness  ratio,  aspect  ratio,  boundary 
conditions,  and  ply  layup  on  the  reactions  of  a  transversely 
loaded  rectangular  composite  plate.  Geometric 
nonlinearities,  i.e.  the  von  Karman  nonlinear  plate 
equations,  and  through  the  thickness  shear  effects  were 
included.  Both  geometrically  linear  and  nonlinear  solutions 
were  done.  This  research  was  done  using  an  existing  finite 
element  code  with  a  foui — noded,  28  DOF  rectangular  element. 

All  plates  were  48  plies  thick  to  be  representative  of 
a  "real  world"  application.  The  ply  layups  were  chosen  to 
be  a  representative  sampling  of  layups  used  in  the  field  and 
had  at  least  SOk  of  the  plies  in  the  0“  direction. 

Material  properties  had  the  largest  effect  when  the 
plate  was  thick.  The  linear  solution  approaches  the 
classical  laminated  plate  theory  solution  as  the  plate  gets 
thin.  The  nonlinear  solution  deflects  much  less  than  the 
linear  solution  because  the  higher  order  terms  significantly 
stiffen  the  plate.  For  a  clamped  plate  the  linear  and 
nonliear  solutions  are  comparible  until  the  plate  is  quite 
thin. 

Aspect  ratio  has  a  large  effect  on  the  plate  for  an 
isotropic  material.  However,  in  the  ply  layups 
investigated,  aspect  ratio  had  much  less  effect  because  half 
the  fibers  were  oriented  in  the  direction  of  the  short 
dimension  of  the  plate. 

Ply  layup  did  not  have  a  great  effect  on  the 
deflections  of  the  plate.  This  occured  because  all  the 
layups  chosen  had  at  least  half  the  plies  oriented  in  the 
same  direction. 
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I.  INTRODUCTION 


Advanced  composite  materials  are  being  used 
increasingly  in  many  engineering  and  civilian  applications, 
from  aerospace  structural  components  to  sports  equipment. 

The  high  stiffness  to  weight  ratio,  coupled  with  the 
flexibility  in  the  selection  of  the  lamination  scheme,  which 
can  be  tailored  to  match  the  design  requirements,  make 
laminated  composite  plates  ideal  structural  components  for 
many  applications.  The  increased  use  of  laminated  plates  in 
various  structures  has  created  considerable  interest  in 
their  analysis. 

In  many  situations,  simplifying  assumptions  can  be  made 
to  help  the  analyst  get  timely  and  accurate  results.  For 
many  structural  components,  only  geometric  nonlinearity  may 
be  important.  For  these  cases,  the  structural  components 
experience  only  small  material  strains  under  load  but  may 
fail  catastrophically  due  to  their  geometric  configuration. 

A  large  class  of  structural  systems  may  be  accurately 
represented  based  on  nonlinear  geometrical,  small  strain, 
and  linear  elastic  material  behavior. 

Many  of  the  nonlinear  displacement  terms  may  be 
considered  negligable  in  a  geometrically  nonlinear  analysis 
depending  on  the  specific  situation.  An  accurate 
load-displacement  characterization  of  a  flat  plate  is  based 
on  the  von  Karman  equations  where  many  nonlinear  rotational 
terms  have  been  discarded  as  negligable.  This 
characterization  gives  an  accurate  representaion  for  plate 
deflections  that  are  many  times  the  plate  thickness.  This 
formulation  is  valid  for  only  moderate  rotations  and  is 
called  an  intermediate  nonlinearity  approach. 


Laminated  plates  also  have  an  additional  complication. 
Three  dimensional  effects  can  become  very  important  because 
of  potentially  large  directional  variations  of  stiffness 
properties  due  to  tailoring.  Also,  the  problem  of 
delamination  consideration,  which  only  now  is  considered  due 
to  low  velocity  impact,  requires  a  higher  order  theory  to 
capture  the  transverse  shear  stress  terms.  Whereas 
classical  two  dimensional  assumptions  may  be  valid  for  an 
identical  plate  structure  consisting  of  isotropic  materials, 
they  may  lead  to  large  errors  for  orthotropic  materials. 
Therefore,  through  the  thickness  effects  should  be 
considered  when  analyzing  composite  materials. 

The  general  class  of  problems  that  include  large 
displacements  and  intermediate  rotations,  small  strains, 
linearly  elastic  material  behavior,  and  parabolic  transverse 
shear  distribution  is  addressed  in  this  thesis. 
Unfortunately,  closed  form  solutions  are  very  limited  and 
cannot  include  all  of  the  desired  features.  Therefore,  the 
finite  element  numerical  approach  is  chosen  as  the  solution 
technique  in  this  research. 

Previous  Work 

Previous  related  work  is  presented  from  two  areas. 
First,  linear  and  nonlinear  techniques  for  analyzing  plate 
structures  are  discussed.  Second,  finite  element  numerical 
solution  techniques  are  discused. 

Plate  Theories.  The  classical  linear  approach,  known 
as  Kirchhoff  flat  plate  theory  C41>  is  based  on  the 
following  assumptions:  1>  the  plate  is  thin,  2)  the 
displacements  and  rotations  are  small,  3)  normals  to  the 
plate  datum  surface  before  deformation  remain  normal  after 
deformation,  and  4)  transverse  normal  stress  is  negligable. 


The  third  assumption  effectively  neglects  transverse 
shear  strains  although  the  transverse  shear  stresses  must  be 
included  in  the  equilibrium  equations.  Since  the  plate  is 
thin,  it  is  assumed  to  be  in  a  state  of  plane  stess,  i.e. 
the  transverse  normal  stress  is  assumed  to  be  negligibly 
small.  Also,  because  the  plate  is  thin,  the  effects  of  the 
normal  tranverse  strains  are  often  neglected  compared  to  the 
effects  of  the  inplane  strains,  i.e.  the  length  of  a  normal 
to  the  datum  surface  is  often  assumed  to  remain  constant 
throughout  deformation.  The  normal  transverse  strain  can  be 
included  in  the  analysis  through  the  constitutive  relations. 
In  deriving  the  equilibrium  equations,  statically  equivalent 
forces  and  moments  acting  on  the  datum  surface  are  defined 
by  integrating  stresses  through  the  thickness. 

As  mentioned  previously,  the  Kirchhoff  flat  plate 
theory  applies  to  a  plate  so  thin  that  all  transverse 
deformation  effects,  i.e.  transverse  stresses  and  strains, 
can  be  neglected.  As  the  plate  becomes  thicker  relative  to 
its  inplane  dimensions,  these  transverse  effects  become  more 
pronounced,  especially  the  transverse  shear  deformations 
CIO, 12, 13).  The  first  theories  that  included  the  tranverse 
shear  deformations  relaxed  the  assumption  on  the  deformed 
normals  to  the  plate  datum  surface.  Now,  the  normal  is 
permitted  to  rotate  so  that  plane  sections  perpendicular  to 
the  datum  surface  remain  plane  but  not  necessarily 
perpendicular.  The  plate  strain  displacement  relations  are 
derived  from  kinematics  and  the  3— D  strain  displacement 
relations.  The  transverse  shear  is  obtained  by  including 
independent  degrees  of  freedom  in  the  kinematics.  The  plate 
is  fully  described  by  the  behavior  of  the  datum  surface  and 
therefore  these  approaches  represent  2-D  theories  C32). 

The  shell  kinematics  developed  by  Bassett,  as  discussed 


in  (8),  express  the  displacefflents  as  an  Infinite  power 
series  in  the  thickness  parameter  or  coordinate. 

Hildebrand,  Reissner,  and  Thompson  CHRT>  C8>  introduced 
truncated  Bassett  kinematics  to  analyze  thin  elastic 
orthotropic  plates  and  the  importance  of  the  transverse 
stresses  and  strains.  HRT  found  that  the  effects  of  the 
second  order  displacement  terms  on  the  transverse  shear 
deformation  were  negligible.  Additionally,  terms  in  the 
transverse  displacement  that  resulted  in  nonzero  transverse 
normal  strains  were  found  to  be  negligable.  Reissner  used 
these  kinematics  to  analyze  plates  C39>.  Mindlin  similarly 
included  rotary  inertia  terms  in  the  dynamic  analysis  of 
plates  C15). 

The  Reissner-Mindlin  CRM)  kinematics  do  not  satisfy  the 
transverse  shear  boundary  conditions  on  the  top  and  bottom 
surfaces  of  the  plate  since  a  constant  shear  angle  through 
the  thickness  is  assumed,  i.e.  plane  sections  remain  plane. 
Because  of  this,  the  theories  based  on  these  kinematics 
usually  require  shear  correction  factors  to  satisfy 
equilibrium. 

Levinson  C14),  Murthy  C16),  and  Reddy  C32>  have 
developed  theories  that  include  cubic  terms  in  the  inplane 
displacement  kinematics  which  produce  zero  transverse  shear 
stresses  at  the  top  and  bottom  surfaces  of  the  plate. 
Satisfying  zero  transverse  shear  stresses  on  the  top  and 
bottom  surfaces  of  the  plate  gives  a  parabolic  shear  strain 
distribution  through  the  thickness,  thus  agreeing  more 
closely  with  linear  elasticity.  The  number  of  variables  in 
the  kinematics  is  equal  to  that  in  the  RM  theory,  but  shear 
correction  factors  are  not  needed. 

The  preceding  discussion  applied  primarily  to  plates 
made  of  isotropic  materials.  HRT  C8)  were  the  first  to 


apply  plate  equal tions  to  an  orthotropic  material.  Reddy 
C32>  discusses  anisotropic  laminated  plate  structures 
including;  various  transverse  shear  deformation  treatments  in 
his  text. 

The  simpifying  assumption  of  laminated  anisotropy  is 
often  used  to  apply  a  2-0  theory  to  plates  made  of  layers  of 
composite  materials  C32>,  In  this  approach,  the  individual 
properties  of  the  composite  elements,  the  fibers  and  the 
matrix,  are  "smeared'*  and  thus  each  lamina  is  treated  as  an 
orthotropic  material.  Also,  laminated  anisotropy  assumes 
perfect  bonding  between  layers,  i.e.  the  interply  adhesive 
has  infinitesimal  thickness  but  infinite  stiffness.  This 
approach  leads  to  classical  laminated  plate  theory  CGLPT>. 
The  references  by  Jones  Cll)  and  Ashton  and  Whitney  c:2>  are 
thorough  presentations  of  CLPT.  GLPT  relies  on  the 
Kirchhoff  assumptions  on  the  datum  surface  normals;  however 
both  references  point  out  that  the  transverse  shear 
deformation  is  more  significant  in  laminated  anisotropic 
structures  over  similar  isotropic  constructions. 

Yang,  Norris,  and  Stavsky  C48)  generalized  the  RM 
theory  to  laminated  anisotropic  plates.  Whitney  and  Pagano 
<47)  were  the  first  to  apply  it  to  composite  plates.  Thick 
composite  plate  closed  form  solutions  were  developed  by 
Reddy  and  Ghao  C35^. 

Reddy  C32,38)  extended  the  cubic  kinematical  approach 
to  analyze  laminated  anisotropic  plates  and  he  and  Soldatos 
applied  them  to  solve  several  linear  static  and  buckling 
problems  C26,33,38). 

Pagano  C20 , 21 , 22 , 231  and  Srinivas  and  Rao  C39) 
developed  some  exact  solutions  of  3-D  elasticity  equations 
governing  composite  plates  that  have  been  used  to  validate 
these  theories.  They  conclude  that  GLPT  gives  fairly  good 
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approximations  for  both  the  displacements  and  stresses  if 
the  plate  is  thin.  Thinness,  as  defined  for  layered 
composite  plates,  not  only  considers  length  to  thickness 
ratios  but  also  the  degree  of  anisotropy.  Transverse 
stresses  are  calculated  from  the  eqi librium  equations  and 
the  CLPT  inplane  stresses.  Higher  order  shear  theories  do 
not  give  much  better  transverse  stress  results  but 
displacements  show  a  marked  improvement  over  CLPT  for  the 
thicker  plates.  Transverse  stresses  are  calculated  best 
from  equilibrium  instead  of  from  the  constitutive  relations 
Cll). 

Much  of  the  previously  mentioned  research  in  the 
analysis  of  composite  plates  is  limited  to  geometrically 
linear  problems.  Reddy  C31>  concludes  this  may  have  been 
due  to  the  complexity  of  the  nonlinear  partial  differential 
equation  associated  with  the  large  deflection  theory  of 
composite  plates.  Approximate  solutions  to  the  large 
deflection  theory  Cin  von  (Carman's  sense>  of  laminated 
composite  plates  were  attempted  by  Whitney  and  Leissa  C46!>, 
Bennett  C4),  Bert  <5),  Zaghloul  and  Kennedy  C49>,  and  Noor 
and  Hartley  C17>.  Zaghloul  and  Kennedy  C49>  used  a  finite 
difference  successive  iterative  technique  in  their  analysis. 
In  all  of  these  studies,  with  the  exception  of  C17>,  the 
transverse  shear  effects  were  neglected.  The  finite  element 
employed  by  Noor  and  Hartley  C17>  includes  the  effect  of 
transverse  shear  strains;  however,  it  is  algebraically 
complex  and  involves  a  large  degree  of  freedom  per  element 
and  thus  one  can  preclude  the  use  of  such  elements  in  the 
nonlinear  analysis  of  composite  plates.  The  finite  element 
developed  by  Reddy  C30)  was  extended  to  nonlinear  bending  of 
composite  plates  C36,37). 

Finite  Element  Plate  Applications.  General  3-D 
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elements  could  be  used  to  model  plate  structures. 
Unfortunately,  since  the  plate  is  thin,  numerical  ill 
conditioning  results.  Also,  using  3-D  elements  would 
involve  many  degrees  of  freedom  CDOF>  that  may  not  be 
necessary.  Consequently,  elements  that  are  specifically 
designed  for  plate  structures  are  developed. 

The  most  popular  approach  is  to  degenerate  the  3-D 
displacement  relations  to  form  the  element.  This  first 
appeared  in  the  form  of  the  Ahmed  element  Cl>.  In  elements 
of  this  type,  transverse  displacements  and  rotations  are 
treated  independently  and  therefore,  are  well  suited  to 
Reissner-Mindlin  shear  deformation  theories.  Also,  normal 
stress  is  assumed  to  be  zero  consistent  with  most  plate 
theories.  However,  it  develops  serious  shear  "locking" 
problems  as  the  plate  becomes  thinner.  Consequently,  this 
method  is  generally  used  with  some  means  to  remedy  the 
locking  phenomenon.  Zienkiewicz,  Taylor,  and  Too  C50> 
introduced  reduced  or  selective  integration  to  alleviate 
these  shear  locking  problems. 

Application  of  material  anisotropy  in  finite  elements 
was  done  by  Pryor  and  Barker  C27>  who  developed  a  linear 
flat  plate  element  based  on  laminated  anisotropic  plate 
theory  including  RN  shear  deformations.  They  suggest  an 
approach  where  each  layer  of  the  laminate  has  rotational 
degrees  of  freedom  which  allows  the  satisfaction  of 
transverse  stress  continuity  at  each  laminate  boundary. 

This  idea  has  been  applied  using  linear  theory  and  the 
Ahmad  element  by  several  investigators  C9).  Recently, 
Palazotto  and  Witt  C25)  extended  the  approach  to  a 
geometrically  nonlinear  shell  formulation  and  applied  it  to 
flat  plates.  Hinrichsen  and  Palazotto  C9)  used  a  cubic 
spline  function  to  represent  the  transverse  displacements  of 
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a  flat  plate  and  found  that  while  ^ivin^  a  hig^her  order 
approximation  to  displacements,  the  equations  were 
simplified  over  the  work  of  C25>  which  essentially  used  a 
quadratic  spline. 

Shear  locking  may  be  reduced  in  finite  element 
formulations  by  using  higher  order  kinematics  to  represent 
the  displacements  through  the  thickness.  Putcha  and  Reddy 
C28^  used  cubic  kinematics  to  formulate  a  mixed  element  for 
nonlinear  anisotropic  plate  analysis.  For  many  of  the  cases 
examined,  the  element  did  not  lock  as  its  thickness  was 
decreased  even  using  exact  integration. 

The  simultaneous  nonlinear  equations  that  result  from 
the  finite  element  discretizations  are  typically  solved 
using  direct  iteration  or  Newton-Raphson  method.  Thorough 
overviews  of  solution  techniques  can  be  found  in  Stricklin 
and  Haisler  <43>,  Bathe  and  Cimento  C3),  Riks  C40>,  and 
Vaszczszyn  C4S).  The  direct  iteration  is  generally  not  used 
because  it  often  diverges,  and  the  Newton-Raphson,  without 
modif ications,  is  not  used  because  it  is  inefficient  and 
sometimes  unreliable.  A  form  of  the  Newton-Raphson  method 
that  is  guaranteed  to  converge  increments  many  small  load 
steps  instead  of  a  single  increment  to  reach  a  target  load. 
One  way  to  make  this  appi'oach  more  efficient,  is  to  not 
update  the  tangent  stiffness  matrix  with  each  iteration. 

This  form  is  called  the  modified  Newton-Raphson  method. 

As  can  be  seen,  much  work  has  been  done  in  the  analysis 
of  composite  flat  plate  structures.  Unfortunaly,  not  as 
much  attention  has  been  placed  on  studying  their  nonlinear 
response.  The  research  effort  here  will  focus  on 
geometrically  nonlinear  composite  plate  behavior  including 
the  very  important  influence  of  transverse  shear 
deformation. 


The  research  presented  here  used  a  finite  element  code 
developed  by  Dennis  C7>.  This  code  was  designed  to  analyze 
a  fully  nonlinear  cylindrical  shell.  However,  since  a  plate 
can  be  modeled  as  a  cylindrical  shell  of  infinite  radius, 
this  code  can  be  used  for  flat  plate  analysis.  The  code 
includes  through  the  thickness  shear  effects  and  the  von 
iCarman  geometric  nonlinearity  equations  as  applied  to 
composite  materials. 


II.  THEORY 


Reference  C7)  contains  a  very  complete  development  of 
the  equations  for  a  finite  element  code  for  a  fully 
nonlinear  cylindrical  shell.  Vi thin  this  code  there  is  also 
the  capability  to  analyze  both  the  linear  and  the  von  Karman 
plate  problems  studied  in  this  research.  Both  types  of 
analysis  include  the  through  the  thickness  shear  effects. 
This  causes  the  linear  solution  to  be  different  than  the 
classical  plate  theory  solution.  This  also  causes  the  von 
ICarman  solution  obtained  to  be  different  than  the 
"classical”  von  Karman  plate  solution.  The  highlights  of 
the  derivation  are  repeated  here  and,  whenever  possible,  the 
derivation  is  related  to  the  von  ICarman  plate  problem.  The 
linear  plate  problem  is  a  simplification  of  the  von  ICarman 
plate  problem. 

Constitutive  Relationships 

Ve  begin  by  assuming  a  linear  constitutive  relationship 
between  stresses,  a  ,  and  strains,  c  ,  of  the  form 

«.  j  '  j 


>•  J 


J kl 


i, j,k,l  =  1,2,3 


C2.  1> 


where  a^^  are  constants  of  the  elasticity  tensor. 

In  this  research,  we  deal  with  materials  that  are 
transversely  isotropic  with  respect  to  planes  parallel  to 
the  2-3  plane.  In  this  case  Eqn  C2.1>  becomes 
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where  contracted  notation  is  introduced 


^l“^ll> 

o  *o 

2  22’ 

°'3*^33’ 

o  =a  .  o  ™o  o  “O 

4  23'  S  13’  6  12 

e  *e 

2  22  ’ 

^3*^33’ 

e_»2e^^,  e_*2e., 

4^  23*  5  13'  ©  12 

C2. 3a) 
C2. 3b) 


Transverse  isotropy  assumes  E  =E  and  and  the 

2  3  X  2  1  3 

^  written  in  terms  of  engineering  constants  become 
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mE*  23  1221  ^  ^  n  •• 

‘  ’  '^23'  ^B5  ^13'  '^60  1: 


■'23  2  A 
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By  assuming  an  approximate  case  of  plane  stress,  i.e.  a  =*0, 

3 

Eqn  C2.2)  becomes 
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For  transverse  isotropy  it  can  be  shown  that 


For  a  plate  constructed  of  layers  of  transversely 
^  orthotropic  material,  the  fibers  are  generally  oriented  at 

some  angle  away  from  the  geometric  axis  of  the  plate. 
Therefore,  the  constitutive  relations  of  Eqn  C2.5>  must  be 
transformed  into  the  plate  coordinates  resulting  in  Eqns 
•  C2.7)  and  C2.8). 
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C2.8) 


Where  the  ??  Ci,j»l,2,6)  and  CF  Cm,n=4,5>  are  elements  of 

V  J  TWT% 

symmetric  arrays  of  transformed  stiffnesses  for  the  k*^  ^  ply- 


Kinematic  Relations 

We  begin  the  derivation  of  the  kinematic  relationships 
by  assuming  a  truncated  power  series  inplane  displacement 
function.  This  will  lead  to  the  desired  parabolic 
transverse  shear  strains.  This  displacement  function  was 
%  developed  in  C28>  and  is  shown  in  Eqn  C2.9)  where  u“,  , 

,  and  are  functions  of  the  inplane  coordinates  only: 


u^  -  u*  + 


C2. 9) 


where  u^  Cl  =1,2)  are  total  inplane  displacements, 
ujis  the  translation  displacement  of  the  datum 
surface, 

<  is  the  through  the  thickness  coordinate, 

is  rotational  displacements  due  to  bending,  and 
B  ,<b  are  unknown  OOF  to  be  determined. 

We  then  assume  that  the  transverse  shear  strains, 
c^=Cg=0  on  the  surface.  For  an  orthotropic  material,  this 
implies  that  the  transverse  shear  stresses,  o^=a  =0  on  the 
surface.  For  example,  for  a  flat  plate  of  thickness  h  with 
^^Ci=l,2>  representing  the  inplane  coordinates 
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Where  the  subscript  ,a  represents  partial  differentiation  of 
the  variable  with  to  respect  to  a.  For  Eqn  <2.10)  to  held 
true  at  the  top  and  bottom  of  the  plate,  Also, 

based  on  <8)  we  neglect  term  in  u  as  small  compared  to 

u“.  Now  solve  Eqn  <2.10)  for 


+ 


] 


<|*±hy'2 


0 


3<h/^)*0j  +  Ug  ^  *  0 

t  -  -  5^.  [  *  n  ] 


<2. 10a) 


This  leads  to  parabolic  transverse  shear  strain  shown  in  Eqn 

<2. 11): 


<2.  11) 


a  similar  approach  can  be  use^  on  e  .  This  leads  to  the 

kinematic  equations  used  in  this  research  as  shown  in  Eqn 

<2. 12); 
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Uj  -  V  *  ^  [  ("a  *  *■  2  ] 


(2. 12a> 

C2. 12b) 
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3  *  * 


C2. 12c) 


The  von  Karman  plate  strain  displacement  relations  are  given 
for  the  inplane  strains  by  Eqn  C2.13): 
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These  come  from  the  Greenes  strain  tensor  C24)  shown  in  Egn 
<2.14): 


2e 

<■  J 


where. 


e  *  Greens  strain  tensor 

>■  j 


2c  »  u  +  u 

J  »■  >  J  J  ,  v 


2o  ■  u  -  u 

>■  J  »•  >  J  J  ,  <■ 


<2.  14) 


<2. 14a) 

<2. 14b) 


Then,  assume  the  following:  1)  strains  are  small,  so  c,  c, 
terms  can  be  neglected,  2)  rotations  relative  to  the  x  and  y 
axis  are  moderately  small,  this  allows  c.  u.  and  a  c 

KtKj  JCLiCJ 

terms  to  be  neglected,  and  3)  rotations  relative  to  <  are 
negligable,  this  allows  o__o__  term  to  be  neglected.  If  one 
Substitutes  Eqn  <2.12)  into  Eqn  <2.13)  and  sets  e  -c  ,  the 
following  is  obtained 
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"''2  ^<^2.2  *  ^  *‘^*•22 


^2,2^  ?"*2 
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"^6  “  “'2  ■"  '''1  -"<<^^1,2  *" 

*  <®k<:2w,^^  +  ,>  +  w,^w. 


where  k 


3h^ 


Notice  that  the  nonlinear  terms  Ci.e.  i-wf  — w?  and  w,  w,  ) 
are  functions  only  of  the  inplane  coordinate  and  not  the 
transverse  coordinate 

If  Eqn  C2.13)  is  rewritten  in  a  more  convienent  form 
for  matrix  manipulations  to  be  conducted  shortly,  the 
following'  is  obtained  where  p  is  handled  using  repeating 
subscript  rule 


where 


where 


IP 


*11  =  ^1,1 

*13  “  kCw,^^  + 

*lp<p-2,4,5,6,7>  .  0 


2p 


e  *  ''»2  7  *'2 

*21  “  ^2,2 

*23  *  *‘^•'22  ■*'  ^2,2^ 


C2. 16> 

C2. 16a> 
C2. 16b) 
C2. 16c) 


C2. 17) 

<2. 17a) 
C2. 17b) 
C2. 17c) 
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*2p<^P-2>4,S,6,7>  -  0 

^6  •  ^6  *  <'’*6p 

where 

“  “»2  *  ^»1  "'1*'2 

*61  “  ^,2  *  'f'2,1 

*63  ”  kC2w,^2  *  *^±^2  *  ^2,1^ 
*6pCp-2,4.5,6.7>  .  0 

Also,  it  can  be  readily  seen  from  Eqn  C2.11>  that 

where 

Eg  -  »,!  +  Vi 
*32  “  3k<t»,j  •*■ 
*gpCp=l,3,4,5,6,7>  -  0 

and 

U  -  <  *  <'’*ip 

where 

®4  *  *'2  ^2 

*42  •  3k<w,2  +  V'2> 

w^pCp=l,3,4,5,6,7)  =  0 


Plate  Potential  Energy 

The  strain  energy  terms  are  given  by 


where 


n  =  U  -  V 

p 


C2. 18> 

(2. 18a) 
C2. 18b) 
C2. 18c) 

C2. 19) 

C2. 19a) 
C2. 19b) 

C2. 20) 

C2. 20a) 
C2. 20b) 


C2.  21) 
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U  ■  Internal  strain  energy 
V  ■  work  done  by  external  forces 

we  divide  the  internal  strain  energy  into  two  parts  Ci.e. 


where  represents  the  inplane  terms  and 
represents  the  through  the  thickness  terms.  Cqns  C2.7>, 
<2.8>,  C2.16),  C2.17>,  <2.18>,  C2.19>,  C2.20>,  and  C2.21) 
give  the  internal  strain  energy  of  the  plate 


"  I  J  J  [  ^11  *  ^22^^2  ^^*2p^^ 

•'h 

+  2  Csl  + 

2  ^26  ^^2  ^%r>  ]  C2.22a> 


U2  =  |-  J  J  [  ^744  CC4  +  <^*42^^  *  ^^*52^^ 

+  2  1743  Ce"  +  <^*42><eg  +  t;%2^  ]  C2.22b> 

where, 

P,r  *  1,2,3 

^  generally  vary  as  a  function  of  <  since  a 
laminate  is  constructed  of  plies  with  different  fiber 
orientations,  and 

ft  represents  the  plate  middle,  or  datum,  surface 

area. 

Also,  if  one  combines  Eqns  C2.16>,  C2.17),  and  C2.18),  the 
following  is  obtained 
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C2.23b) 
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Ve  can  now  represent  the  first  part  of  the  internal  strain 
energy  more  convieniently  thus 


i  r  r  CQsy'^e  d«;dfl 

•^0  -^h 


C2. 24> 


Eqn  C2.241  can  be  rewritten  in  terms  of  the  plate  mid 
surface  where  the  <  dependence  has  been  integrated  out  by 
defining  a  series  of  elasticity  arrays 


where. 
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where  1,J  ■  1,2,6  and  p,r  *  1,2,3  Crepeatin^  subscript 
rule  again  applies>  and 


[A  ,B  ,D  ,E  ,F  ,G  ] 

'-'•J  V-J  l-J  >-J  v-J  '•J-' 


C2. 25d) 


Substituting  Eqns  C2.16^,  C2.17>,  and  C2.18>  into  Eqns 
C2.25a),  C2.25b),  and  C2.25c>  results  in  the  following: 
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C2. 26> 

“2 

-  2  Cj  ^*il®ij  *i3^ij^ 

C2.27) 
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“  *jl*il*^ij  *  2*jl*i3*^ij 

C2. 28> 

The  von  Karman  plate  relations  only  retain  the 
nonlinear  displacement  terms  in  the  transverse  displacement, 
w.  Note  that  only  cubic  and  quartic  terms  will  result  in 
nonlinear  terms  in  the  equilibrium  equations  C<5np*0>.  The 
only  cubic  and  quartic  terms  in  w  come  from  C2.26>  and 
C2.27>.  Also,  for  symmetrically  arranged  laminates  Ci.e. 
the  case  in  this  research^  ^  i“  ® 

A  J  A  vl 

Therefore, 


U.  -  if  e’.  cTA.. 


dfl 


C2. 29> 


and  with  a  similar  treatment  as  ,  U2  becomes 
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C2.  30:> 


U~  »  i  r  CC*C*A„„  +  ■*■ 

2  ^  0  I**  ^  ^  ***^  ***^  ***^ 

where 

C^n'  D„„,  F^3  -  /^[l,  <2_  jj. 

Now  we  have  all  the  pieces  in  place  for  the  "modified”  von 
ICarman  approach.  Modified  in  that  it  includes  the  through 
the  thickness  shear  terms. 

Finite  Element  Solution 

The  first  part  of  the  derivation  is  independent  of  the 
particular  element  chosen.  The  second  part  concerns  the 
specific  element  used  in  this  research.  We  begin  the 
element  independent  portion  of  the  derivation  by  developing 
an  expression  for  the  strain  energy  of  the  form 


r  N.  N  , 

"p  *  J  [  •=  *  -3^  *  ^ 

where, 

q  is  a  column  array  of  nodal  displacements, 

R  is  a  column  array  of  nodal  loads, 

K  is  an  array  of  constant  stiffness  coefficients, 
N^is  an  array  of  stiffness  coefficients  that  are 
linear  in  displacement,  and 

N^  is  an  array  of  stiffness  coefficients  that  are 
quadratic  in  displacement. 

Taking  the  first  variation  if  Eqn  C2.31>  gives 


*  <5q'*'  FCq)  -  0 


C2.  32> 
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FCq)  represents  the  equilibrium  equations  For  an  arbitrary 
and  independent  <5q 

FCq)  -  0  C2.33> 

FCq)  -  4  K  +  +  --2.  J.  q  -  R  -  0  C2.34> 

Eqn  C2.33)  can  be  expanded  in  a  Taylor  series  giving 

aF 

FCq  +  Aq)  “  FCq)  +  ^  Aq  +  «  o  C2.35) 

Now  neglect  higher  order  terms  because  Aq  is  small  and 
rearranging  gives 

—  Aq  -  -FCq)  C2.36) 

aq 

[  K  +  N^  +  Njj  3  Aq  =  -FCq) 

Aq  -  -FCq)  C2.3r) 

where,  “  tangent  stiffness  matrix 

K  =  K  +  N,  +  C2.38) 

T  12 

Eqn  C2.37)  is  solved  in  an  iterative  manner  by  a 
Newton-Raphson  technique.  The  current  values  of  q  are 
inserted  in  N^and  N^  on  both  sides  of  Eqn  C2.37)  and  the 
ensuing  linear  equations  are  solved  for  Aq.  Aq  is  then 
added  to  q  and  the  new  q  is  inserted  into  N^and  and 

so  on.  This  is  continued  until  the  RHS  becomes  arbitraily 
small . 

The  proceeding  derivation  assumed  a  definition  of  K, 

N^,  and  N^,  such  that  these  terms  repeated  themselves  in  the 
first  variation  Ci.e.  Eqn  C2.31)  to  Eqn  C2.32))  and  aqain  in 
the  Taylor  series  expansion  Ci.e.  Eqn  C2.33)  to  Eqn  C2.35)). 
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Rajasekaran  and  Murray  C29)  and  Dennis  <7>  have  developed  a 
^  procedure  which  insures  this  repetition.  The  overall 

technique  is  shown  here  for  a  "modified"  von  Karman  flat 
plate.  For  a  more  detailed  development  see  the  references 
listed  above. 

#  Ve  be^in  by  dividing  each  strain  component  of  Eqns 

C2.1d>,  C2.17>,  and  C2.18)  into  linear  and  nonlinear  parts: 


®i  ■  •“! 


lT  d  +  i.  d'*’  H,  d 


•  1 


ip  pi 


lT  d  +  i  d"^  d 


P  1 


C2. 39) 


where  “  0 , 1 , 2, 3 , 4 , 3 , 6 , 7)  are  column  arrays, 

•H^  are  symmetric  arrays, 
i  *  1,2,6  and  p  =■  1  to  7,  and 
d”*'  is  the  displacement  gradient  vector. 

d”^  =  {  u  u,^  u,2  V  v,j  v,2  w  w,^  w,2 

’'’22  ’’’12  n  ^1,1  ^1,2  <^2,1  ^"2,2  }  ^2.40) 


Each  .L^are  18x1  column  arrays  and  each  {H^is  an  18x18 
array.  From  Eqns  C2.16),  C2.17),  and  C2.18)  it  can  be  seen 
that 


oW^2)  =  1 


L,C6)  =  1 

O  2 


„L.C3)  =  „L^C5)  =  1 
o  e  0  6 


jLjC14)  »  1 


L,C18)  -  1 
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1  O  X  o 


C2. 41a) 
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3L/10>  -  3L/14>  -  k 

3L3CII)  -  3L3CI8)  -  k 

3L3C12>  -  2k,  3L3C15>  -  3L3Cir>  -  k 

H  :  C8,8>  -  1 

01  ' 

C9,9>  «  1  C2.41b) 

02  ' 

3H3:  <8,9>  -  2 


where  the  number  in  the  parenthesis  denotes  the  position  of 
the  given  value  in  the  particular  array.  All  other 
positions  are  zero. 

Now  substituting  Eqn  C2.39)  into  Eqn  C2.29>,  which  is 
repeated  here,  one  obtains  Crecall  these  expressions  are  for 
inplane  strain  energy  using  strain  at  the  datum  surface) 


gives 


where 


U.  =  1.  f  c\  cT  A.  .  dfl  C2.29) 

1  2  J  i  ij 


-  r 
*  -^0 

d’’ 

C2. 42) 

a  - 

O^v 

0  j 

C2. 43a) 

*ij 

oL. 

oH 

0  j 

C2. 43b) 

0«v 

d  d’*'  H 

0  J 

C2. 43c) 

%  and  the  repeating  subcript  rule  again  applies. 

Unfortunately,  not  all  these  terms  show  the  required 
repetition  when  the  variation  and  derivative  are  taken. 
That  is,  new  terms  will  result.  In  order  to  get  required 
repititions,  the  following  substitutions  are  made: 

In  Eqn  C2.43b)  let 
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L  <!■*■  H  »  i  Tl  +  dT^L  H  +  H  d  L  1  C2.44) 

Oi  Oj  3|_Oi  Oj  Ovoj  Oi  OjJ 

and  in  Eqn  C2.43c>  let 

d  d'*'„H  -  1  Lh  d  d'*'„H  +  i  d'*‘„H  d  + 

0\.  Oj  3^0  1  Oj  2  Ot  Oj 

d  d'*'^H  +  i  d***  H  d  H  1 

Ot  Oj  2  oi  ojJ 

With  these  substitutions,  Eqn  C2.42>  can  be  rewritten 


C2.45) 


N 


1  + 

3  ”6 


A] 


d  dO 


C2. 46) 


The  tranverse  shear  energy  terms  are  handled  in  a 
similar  manner.  These  are  much  simpler  because  the 
tranverse  shear  terms  are  linear  in  displacement  and 
therefore  all  H  terms  are  zero.  Therefore,  let 


j  m 


-  s’*"  d 
m  O  m 


d 

^m  2  m 


m  ®  4 , 5 


C2. 47) 


From  Eqns  C2.19)  and  C2.20)  it  can  be  seen  that 


oS^^9) 

-  oVl^> 

ae 

1 

-  „S^C13) 

1 

• 

2S^C9) 

=  ,S^C16) 

s 

3k 

2V8> 

=  ,S^<13) 

- 

3k 

C2. 48) 


Because  these  strains  are  assumed  to  be  linear,  strain 
energy  terms  brought  about  by  through  the  thickness  shears 
are  only  given  for  H.  and  are  zero.  Thus,  by 
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substituting  Eqn  C2.47>  into  Eqn  C2.30>  one  obtains  the 
internal  energy  function 


“^•4-  r  A  S  S’’+2D  ,S'*'+F  S  ,S‘*'ld  dft 

Z  2  I  I  mn  O  m  O  n  mn  O  m  2  n  mn  2  rn  2  n  I 
O 

“  C2.49) 

The  repeating  forms  for  the  first  variation  and  Taylor 
Series  expansion  which  are  added  to  1C  are  formed  by  making 
the  substitution  below: 


s  s'*” 

p  m  r  n 


1 

2 


s 


r  1 


s 

r  t 


'T] 


p,r  »  0,2 


C2.  50> 


The  last  step  in  the  element  independent  portion  of  the 
derivation  is  the  discretation  of  the  domain.  This  is  done 
by  approximating  the  continuum  displacements  by 
interpolation,  or  shape  functions.  So  let 


u  =  yf  q 


C2. 51> 


where, 

u  *  vector  of  continuum  displacements 
/f  *  array  of  shape  functions 
q  =>  nodal  values  of  displacements 
Then,  the  displacement  qradient  vector  d  is  approximated 
from  Eqn  C2.51)  by 


d  ■  ©  q 


C2. 52) 


where  ©  is  the  array  of  shape  functions  and  their 
derivatives.  ©  is  defined  for  a  specific  element  and  will 
be  derived  subsequently.  Substituting  Eqn  C2.43)  into  Eqn 
C2.46)  and  including  the  tranverse  shear  terms  of  Eqn  C2.49) 
in  K,  the  following  potential  energy  function  is  obtained 


a'*’  r  ^2  1  T 

“p  -  I  [  ♦  -3^  -  ]  1  - 


C2. 53) 


In  Eqn  C2.53) 


K  -  f  1C  2>  dO 

•'o 


N  «  r  N  ®  dft 
1  1 


C2. 54) 


N"  r  ©■*■  N,  2>  dO 
2  2 

Derivation  of  g  for  a  28  DOF  element 

Now  we  begin  the  second  portion  of  the  finite  element 
derivation.  This  part  is  element  dependent.  The  continuum 
displacements  are  approximated  by  the  nodal  displacements 
and  shape  functions.  These  shape  functions  are  chosen  such 
that  the  following  requirements  are  met; 

Cl)  continuous  displacements  within  an  element , 

C2)  the  element  must  be  able  to  represent  constant  strain, 
C3)  rigid  body  modes  are  present  in  the  assumed 
displacements, 

C4)  compatibility  exists  between  elements,  and 

C5)  the  element  should  not  have  a  preferred  direction,  this 

is  desired,  not  required. 

From  our  displacement  function,  Eqn  C2.12),  it  can  be 
seen  that  seven  degrees  of  freedom  CDOF)  are  present  at  each 
node.  These  are  u,  v,  w,  w,^,  w,2»  ,  and  V'2-  Continuity 

of  displacement  only  is  required  for  u,  v,  ,  and  v'2- 
Therefore,  Lagrangian  shape  functions  may  be  used  since  they 
gaurantee  G^continuity .  Continuity  of  displacement  and 
their  derivatives  is  required  for  w  Ci.e.  w,  *’2^ 
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Hermltian  shape  functions  can  be  constucted  that  will 
provide  c‘  continuity. 

Dennis  C7>  begins  the  element  definition  by  assumini; 
the  following  12  term  quartic  for  the  transverse 
displacement  w,  of  a  plate  with  four  corner  nodes.  See 
figure  2.1. 

wCx,y)  “ 

+  a^x^y  +  a^xy^  +  ^lo^^  *  C2.55) 


X 


FIGURE  2.1.  Rectangular  Plate  Element  with  Four  Nodes 

Now,  Eqn  C2.55>  can  be  rewritten  in  the  form  of  Eqn  C2.51) 


for  an  element  as  shown  below. 


1x3 


wCx, y> 


[  *3  ] 

1x12 


1x3 

. 

^1 

^3 

12x1 


C2. 56) 


Cnote:  Cn  x  m)  is  the  matrix  order  indicator) 
where,  for  the  k*"  ^  node. 


<  w  C2.56a) 

Hermitian  shape  functions  shown  below 


f  •> 

< 

k  1 

at, 

k  2 

k. 

ici +<  ^  <  )  <  1 +17,^  T?)  C  2+C  ^  < +r7  ^  17-t;  ^ -r)^  ) 


C 1 +<  ^  <)  C  rr-l )  C 1 +77k  »?>  * 


C2. 56b) 


node  has  the  natural  coordinates  CC.^,rf^'>  and 


a-'k'*  ''k'^'  '‘'k' 

where  the  k' 
i:;=x/a  and  rj^y^h.  See  figure  2.2. 

For  the  remaining  nodal  DOF  u,  v,  and  14^2  only  C" 

continuity  is  required  so  we  will  use  Lagrangian  shape 
functions.  For  these  DOF,  we  approximate  Eqn  C2.51)  in  the 
form  shown  below: 


'  U  ' 

r 

0 

0 

0 

... 

0 

0 

0  ■ 

r  ^1  1 

V 

-  ® 

0 

^1 

0 

0 

•  •  • 

0 

J¥ 

4 

0 

0 

•« 

^z 

1 

0 

0 

^1 

0 

•  •  • 

0 

0 

0 

^3 

1  2  J 

•4x  1 

0 

0 

0 

^1 

#  •  • 

0 

0 

0 

1  J 

^x 1 e  1  Ox  1 

f 

where  for  the  k  node. 


qj  -  <  u  V  v/,  V^2  >, 

are  Lagrangian  shape  functions 

^k  -  7 


C2. 57) 


C2. 57a) 


C2. 57b) 
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FIGURE  2.2.  Rectangular  plate  element  of  Figure  2.1  in 
natural  coordinates,  <=X/^a  and  rj^y^b.  This 
is  the  only  element  used. 


Combining  Eqn  C2.56>  and  Eqn  C2.57>,  the  displacement 
gradient  vector,  Eqn  C2.40>,  may  be  expressed  as  in  Eqn 
C2.58)  below: 


dC<,i7>  -  Dq 

18x1 


■  INj  0  0 

»  •  • 

0  0  ■ 

0  IH  0 

•  •  • 

0  DH  0 

- 

; 

L  0  0  (N^ 

•  •  • 

0  0 

1  8x2  8  2  8x 1 


C2. 58) 


where  for  the  k*"  node. 


<  u  V  w  w,  ^ 


^1  >) 


C2. 58a) 
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and 


■ 

0 

■  *ki 

*k2 

*k3 

^k,C, 

0 

*kl,< 

*k2,< 

*k3,< 

0 

*kl,r7 

*k2,i7 

*k3,i7 

0 

^k 

*ki,<t; 

^k2,« 

*k3 , « 

0 

*kl ,  1717 

^k2,rjrj 

*k3,r7T7 

0 

^k,n 

*kl,77< 

*k2,r7< 

*k3,i7< 

6x2  6x3 


C2. 58b> 


The  displacement  gradient  vector  in  terms  of  elemental 
coordinates  <x,y)  of  Figure  1  is  related  to  the  element  in 
natural  coordinates  C<,j7)  Ci,e.  Eqn  C2.S8)  and  Figure  2>  by 
the  Jacobian  matrix,  J.  A  general  derivation  of  the 
Jacobian  follows.  Assume  a  function  fCx,y>  that  represents 
the  elements  of  d  in  Eqn  C2.40>  where  x,y  are  functions  of 
the  natural  coordinates  Taking  the  partial  derivatives 

of  f  gives 


then. 


C2. 59) 


C2.  60) 


In  addition,  second  derivatives  are  required  of  the  w  DOF, 


• 

for  example. 

from  Eqn  C2. 

59) 

rx 

II 

X  '  K, 

+  Cf, 

X 

• 

^^’xx  * 

f 

’  xy 

y,^)x,^  + 

X 

+ 

^  ^ ’ yy  ^ ’ K  * 

f 

’  xy 

x,^)y,^  + 

f'.y 

C2. 61) 
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Usln$  similar  expressions  for 


r  f  1 

'  XX 

r  X* 

A  >  ^ 

•< 

r 

’yy 

X* 

• 

f 

L  ^ ' xy  J 

x,^x. 

0 

■ 

■ 

o 

— 

^>Kri 

^'<0  . 

- 

2y  > 

2y,^x,^ 


ty»  ^y. 


-1 


r  r 

- 

f , 

'ITT? 

I 

C2.62)  i 


Rearrangini^  Eqn  C2.62>  and  combining  terms  into  a  3x5 
matrix.  A,  gives 


r  r  'I 

f  f  ^ 

^xx 

^’yy 

.  .  [  A  ]  . 

.  ^'xy  , 

f  „ 

1  J 

C2. 63> 


#  Combining  Eqn  <2.60)  and  <2.63)  gives  the  inverse  of  the 

Jacobian  matrix  and  relates  the  coordinate  systems  if 
Figures  1  and  2  by  the  following 


where. 


and 


d<x,y)  »  r  d<<;,>7) 


fi  0 

0  ^2 

0  0 

18x18 


0 

0 

r. 


<2. 64) 


<2. 64a) 
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■  1 

0 

0 

0 

0 

0 

• 

0 

Til 

^12 

0 

0 

0 

0 

r2i 

^22 

0 

0 

0 

Tl  - 

0 

0 

0 

1 

0 

0 

0 

0 

0 

0 

^11 

^12 

• 

0 

0 

0 

0 

^21 

^22 

Ox  6 

■  1 

0 

0 

0 

0 

0 

• 

0 

Til 

ri2 

0 

0 

0 

0 

^21 

^22 

0 

0 

0 

^2  - 

0 

^1 

^^12 

'^13 

^4 

^5 

• 

0 

^21 

^22 

^23 

^24 

^^25 

0 

^31 

^2 

^33 

^34 

^35 

6x  e 


C2. 64b> 


C2. 64c) 


For  the  rectangular  element  assumed  here,  many  terms  are 
zero.  What  we  are  left  with  is  r22“l^b>  A^3=l/'a*, 

A2^*1/^*,  and  Agg^l/^ab.  Where  a  and  b  are  the  dimensions  of 
the  plate  in  the  x  and  y  direction,  respectively.  The 
remaining  f's  and  A’s  are  zero. 

The  28  DOF  element  shown  in  figure  2.3  is  now 
completely  defined  and  is  summarized  below. 


d  »  r  D  q  C2.65) 

=  ©  q 

where  the  order  of  the  appropriate  matrix  can  be  stated  as 


d  =  18  X  1 
r  =  18  X  18 
D  =  18  X  28 
q  =  28  X  1 


CEqn  C2.40>> 
(Eqn  C2.64a)) 
CEqn  C2.58) 
CEqn  C2. 58a> 


©  is  the  matrix  that  is  used  in  Eqn  C2.54)  to  calculate 
the  components  of  the  tangential  stiffness  matrix. 
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FIGURE  2.3.  28  DOF  Rectangular  Plate  Element. 


Solution  Algorithms 

Now  that  all  quantities  have  been  defined,  the  solution 
algorithm  is  presented.  The  integrals  inherent  in  Eqn 
C2.37>  are  calculated  for  an  individual  element.  Each 
element's  contribution  is  added  to  the  global  stiffness 
array  according  to  the  element's  connectivity.  Therefore, 
Eqn  C2.37)  is  represented  below  for  the  entire  plate. 


n 

1 

k  =1 


f  2>'*' 

K  +  N 

1  2>  dQ„ 

No 

L  1 

2 

1  ® 

Aq  - 
k 


C2.  66> 


n 


-  5 

k  =1 


I 

r  K  +  A  ^ 

K 

]  ® 

L  2 

3 

J 

where , 

=  2— D  domain  of  the  element 
e 

n  number  of  elements  in  mesh 


*  global  column  arrays  of  displacement  assembled 
Troffl  element  Aq  and  q. 

K  «  applied  load  array^  or  same  dimension  as  global  q. 

R  can  be  found  using  the  same  shape  functions  as  used 
for  element  displacements  and  can  be  calculated  an  an 
elemental  Cand  then  assembled)  or  global  basis.  A  dicussion 
of  both  methods  can  be  found  in  Cook  (6).  Solutions  to  Eqn 
C2.66)  are  found  iteratively  where  the  RHS  becomes  small, 
i.e.  equilibrium  is  satisfied.  The  integrations  are  carried 
out  using  Gaussian  quadriture  in  natural  coordinates  <  and 
17.  As  an  example,  look  at  the  first  term  in  the  summation 
of  Eqn  C2.66): 


J  K  +  » 


1  1 


J  J  [  K  +  Nj  +  ®  C2. 67> 


-1  -1 


where , 


''i  ^  ”1  ^2  ]  ®  IJI 


|J|  =  determinant  of  the  Jacobian  matrix, 

©,  K,  N^,  |J|  are  evaluated  at  Gauss  integration 


points  Ce;.  5, 


W. ,  *  weighting  factors,  and 

sJ 

m  =  order  of  the  numerical  integration. 


In  a  linear  analysis,  displacements  and  rotations  are 

^  A 

assumed  to  be  small,  and  therefore  and  are  eliminated. 
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Also,  since  the 
solved  directly 
necessary.  So, 


relationships  are  linear,  Eqn  C2.34>  is 
for  q,i.e.  iterative  solutions  are  not 
Eqn  C2.66>  may  be  rewritten 


n 


2 

k  =1 


J  J  ^ 


]  2>  1J|  d^drj 


q  -  R 


Jk 


C2. 68> 


The  solution  of  Eqn  C2.68>  is  accomplished  by  Gaussain 
elimination.  The  resulting  displacements  can  be  used  to 
find  the  strains  and  then  stresses  at  the  Gauss  points  via 
previously  derived  relationships. 

For  a  nonlinear  analysis,  the  linearized 
incremental/iterative  equations  are  solved  by  a 
Newton-Raphson  technigue.  For  this  research  5x5  integration 
was  used  rather  than  the  7x7  called  for  in  the  exact 
solution  for  the  nonlinear  cases  and  4x4  integrations  for 
the  linear  cases.  Dennis  <7>  determined  that  this  gave 
adequate  accuracy  and  accelerated  computer  run  times 
considerably.  In  the  first  iteration  of  the  first 
increment,  the  displacement  gradient  vector  is  calculated  by 
the  linear  technique  described  in  the  preceeding  paragraph. 
N,  and  are  then  calculated  and  iteration  continues  until 
convergence  is  achieved.  A  global  displacement  convergence 
criterion  is  used  C19)  as  shown  below: 


X  100  <  TOL 


C2. 69> 


where. 


^ ,  q^  are  the  elements  of  q  for  the  r*’'^. 
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and  first  iterations  of  a  ^iven  increment,  and 
TOL  is  a  user  defined  convergence  tolerance  percentage 
and  was  used  in  this  research. 
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III.  NUMERICAL  DEVELOPMENT 


In  this  section  the  results  of  this  research  effort  are 
presented.  The  results  are  presented  in  four  major  areas. 
First,  the  results  of  the  code  verification,  second,  the 
specific  parameters  of  the  cases  studied,  third,  the 
generation  of  the  finite  element  mesh,  and  finally,  some 
considerations  reguarding  symmetry  of  the  mesh. 

Code  Verification 

The  finite  element  code  used  in  this  research  was 
developed  by  Dennis  C7)  for  a  fully  nonlinear  cylindrical 
shell.  It  can  also  be  used  for  the  plate  under 
consideration  in  this  thesis.  The  code  was  verified  as  to 
its  accuracy  by  Dennis  C7).  The  verification  presented  here 
is  not  only  to  retest  the  accuracy  of  the  code,  but  also  to 
assure  that  the  input  decks  of  the  program  were  formulated 
properly.  Two  verification  cases  were  studied:  a 
geometrically  linear  analysis  and  a  geometrically  nonlinear 
analysis. 

For  the  geometrically  linear  case  the  orthotropic 
material  studied  was  aragonite  crystals  with  the  following 
properties  given  in  million  psi: 


20. 83, 
3.71, 


E^  -  10.94, 

Q^3  =  6.19,  1.^3  =  0.44,  =  0.23 


C3. 1) 


The  plate  was  rectangular  in  dimension  a  x  b  with  a 
uniformly  distributed  transverse  load.  The  transverse 
deflection  w  in  table  3.1  is  taken  at  the  center  of  the 
plate,  h  is  the  thickness  of  the  plate,  and  is  the  load 
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per  unit  area  of  the  plate.  The  values  for  comparison  come 
from  the  exact  solutions  of  Srinivas  and  Rao  C42>  as 
presented  by  Reddy  C32>. 


Table  3.1.  Verification  of  geometrically  linear 
analysis. 


a/b 

a^h 

Exact'’ 

Present 

7. 142 

387. 23 

387.90 

O.S 

10.0 

1408. 5 

1410.0 

20.0 

21542 

21552 

7. 142 

191.07 

191.71 

1.0 

10.0 

688.75 

689. 83 

20.0 

10443 

10447 

7.142 

39.790 

40.240 

2.0 

10.0 

139.08 

139. 90 

20.0 

2048.7 

2051.8 

=  23.2  X  10®  psi 


•  ^Srinivas  and  Rao 


As  can  be  seen,  there  is  excellent  agreement  between  the 
exact  and  the  present  finite  element  solution. 

For  the  geometrically  nonlinear  case  a  test  done  by 
Dennis  C7)  was  repeated.  Dennis,  in  turn,  was  verifying  the 
results  of  Putcha  and  Reddy  C28).  In  this  case  the  plate 
was  16  inches  square,  1.6  inches  thick,  and  subject  to  a 
uniform  pressure  q^.  The  ply  layup  was  [0/'±45/90]^  with  the 
following  material  properties  Cwhere  the  modulii  are 
expressed  in  million  psil: 

E^-60,  E^-l.S,  .9,  G^3-.75,  C3.2> 
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Both  simple  and  clamped  supports  with  the  boundary 
conditions  as  shown  in  Eqn  C3.3>,  were  investigated: 


Simple 

@x»±a-<^ :  v«w“V'2*0 

@y  ■±a/2 :  u-w-v^  j 


Clamped 

U«W*W, 

v-w-w,  2“V'2*0 


C3. 3> 


The  pressure  and  center  displacement  were  nondimensionalized 
as  in  Eqn  C3.1>: 


q  =■ 


E^h^ 


C3.  4) 


The  plate  was  modeled  using  quarter  plate  symmetry  with  a 
4x4  mesh  Tor  the  simple  case  and  8x8  mesh  for  the  clamped 
case.  Results  are  shown  in  Table  3.3. 


Table3.2.  Verification  of  geometrically  nonlinear  analysis 


Dennis 

Present 

M 

w  CSimple) 

w  CClamped> 

w  CSimple)  w 

CClamped> 

so 

.2759 

.  1395 

.  2759 

.  1395 

100 

.  4873 

.  2593 

.  4873 

.  2595 

130 

.  6473 

.  3373 

.  6473 

.  3573 

200 

.  7753 

.  4377 

.  7753 

.  4377 

250 

.  8825 

.5055 

,  8827 

.  5055 

As  can  be  seen  from  table  3.2  the  results  match  extremely 
well.  Also,  the  comparison  of  these  results  with  Pucha  and 
Reddy  can  be  seen  in  figure  3.1.  Again,  agreement  is 
excellent. 

Parameters  Studied 

In  this  research  two  types  of  material  were  studied. 


FIGURE  3.1.  w  vs  q  for  [0/±45/90]a  Laminated  Plate. 
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alon^  with  two  types  of  analysis,  six  thickness  ratios, 
three  aspect  ratios,  two  boundary  conditions,  and  four  ply 
layups.  All  possible  combinations  of  these  six  parameters 
were  analyzed  for  a  total  of  576  test  cases. 

Materials.  The  two  types  of  materials  were  chosen  to 
have  a  large  difference  in  their  ratio  of  E,  to  This 

ratio  is  of  interest  because  it  is  one  of  the  driving 
parameters  in  determining  how  much  effect  the  transverse 
shear  stresses  have  on  the  response  of  the  structure.  The 
specifics  of  each  material  type  are  shown  in  Eqn  C3.5); 


Material  A: 


1 

'2 


40, 


G 


1  3 


0.5 


G 


23 


0.2, 


^12=  ‘^i3“0-25 


C3. 5a) 


Material  B: 


1  _ 


=  15, 


13 

'2 


=  0.6, 


G 


2^3  _ 
2 


*  0.3, 


1  2 


3-0. 28 


C3.  5b) 


For  both  materials  6^2“®! 3  thickness  was  set  at 

0.005  inches.  For  each  material  E^-20xl0®  psi. 

Analysis.  The  two  types  of  analysis  conducted  were 
geometrically  linear  and  geometrically  nonlinear.  Both 
types  of  analysis  were  discussed  in  the  previous  chapter. 

In  early  verification  tests,  it  was  seen  that  loads  greater 
than  500  psi  would  cause  large  enough  rotations  so  as  to 
violate  the  von  Karman  limitations  of  moderate  rotation,  the 
criterion  being  that  the  rotations  could  not  be  greater  than 
17°.  Therefore,  for  both  types  of  analysis  a  maximum  load 
of  500  psi  was  applied  to  the  plate.  In  the  linear  cases 
the  load  was  applied  in  one  increment.  In  the  nonlinear 
cases  the  load  was  applied  in  five  equal  increments  of  100 
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psi.  The  nonlinear  computer  run  times  were  roughly  20  times 
greater  than  the  linear  computer  run  times. 

Thickness  Ratio.  Thickness  ratio,  S,  is  defined  as  the 

length  of  the  plate  in  the  x-direction,  a,  divided  by  the 

thickness.  Six  thickness  ratios  were  studied.  They  were 
10,  20,  30,  40,  50,  and  60. 

Aspect  ratio.  Aspect  ratio,  AR,  is  defined  as  the 
length  of  the  plate  in  the  x— direction,  a,  divided  by  the 

length  of  the  plate  in  the  y— direction,  b.  Three  aspect 

ratios  were  investigated:  0.5,  1,  and  2. 

Boundary  conditions.  The  two  boundary  conditions 
considered  were  simple  and  clamped.  The  boundary  conditions 
were  chosen  so  as  not  to  restrict  the  inplane  deflections  at 
the  boundary.  This  was  done  to  eliminate  as  much  of  the 
inplane  membrane  stresses  as  possible  so  the  effects  of  the 
transverse  shear  stresses  would  be  more  pronounced.  The 
boundary  conditions  are  shown  in  figure  3.2  where  subscript 
t  denotes  tangential  direction  and  subcript  n  denotes  normal 
direction.  In  both  cases  displacements  in  and  out  of  the 
plane  of  the  page  are  allowed. 

Ply  layups.  The  four  ply  layups  chosen  were: 
to,  /90^]  ,  I0^y60y-601  ,  to, ,/45  /-45^]  ,  and 
to, „/45,/-45,/90, ]  ,  where  the  numerical  subcript  denotes 
the  number  of  plies  in  a  given  orientation  and  the  subscript 
s  denotes  that  the  ply  layup  is  symmetric  about  the 
midsurface.  These  layups  represent  several  of  the  more 
common  types  of  ply  layups  used  today.  The  ply  layups  were 
all  48  plies  thick,  which  is  of  the  same  order  as  the  number 
of  plies  used  in  aircraft  structures  such  as  wing,  tail,  and 
stabilizer  skins.  Additionally,  it  was  desired  to  have  at 
least  50k  of  the  plies  in  the  0 ^  or  x-direction. 
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Figure  3.2.  Boundary  conditions. 

Finite  Element  Mesh  Generation 

Because  the  number  of  plies  was  held  constant  and  the 
ply  thickness  was  the  same  for  both  materials,  the  thickness 
of  the  plate  was  constant  for  all  cases  studied.  All  plates 
were  0.24  inches  thick.  Consequently,  the  thickness 


determined  the  inplane  dimensions  of  the  plate  through  the 
thickness  ratio  and  aspect  ratio  of  interest.  The  dimension 
in  the  x-direction,  a,  was  set  by  the  thickness  ratio 
desired.  For  example,  if  a  thickness  ratio  of  30  was 
needed,  then  the  a  dimension  would  be  7.2  inches.  The 
dimension  in  the  y-direction,  b,  was  then  set  by  the  aspect 
ratio.  Continuing  the  above  example,  if  a  was  7.2  inches 
and  an  aspect  ratio  of  2  was  desired,  then  b  would  be  3.6 
inches.  The  values  of  a  and  b  for  all  thickness  and  aspect 
ratios  is  shown  in  table  3.3. 


Table  3.3.  Dimensions  of  plate  for  all 
thickness  and  aspect  ratios. 


s 

a 

bCAR®. 5> 

bCAR=l> 

bCAR»2) 

10 

2.  4 

4.8 

2.  4 

1.2 

20 

4.8 

9.  6 

4.8 

2.  4 

30 

7.2 

14.4 

7.2 

6.6 

40 

9.6 

19.  2 

9.  6 

4.8 

30 

12.  0 

24.0 

12.0 

6.0 

60 

14.4 

28.  8 

14.4 

7.2 

To  illustrate  the  dimensions  of  the  plate  for  any  given 
value  of  S,  the  relative  sizes  of  the  plate  for  the 
different  aspect  ratios  are  shown  in  figure  3.3. 


FIGURE  3.3.  Relative  plate  planforms  Tor  any  value  of  S. 


For  all  cases  studied,  an  8x8  uniform  mesh  was  used.  Dennis 
C7)  showed  that  this  gives  good  accuracy  with  efficient 
computer  run  times.  Also,  quarter  plate  symmetry  was  used 
in  all  cases.  Quarter  plate  symmetry  allows  for  a  finer 
mesh  over  a  smaller  portion  of  the  plate  and  thus  increases 
accuracy  while  decreasing  computer  run  time.  Dennis  C7) 
also  ran  a  convergence  test  which  verified  the  accuracy  of 
the  8x8  quarter  plate  mesh. 

Quarter  Plate  Symmetry  Considerations 

Early  in  the  research  effort  the  question  arose  as  to 
whether  quarter  plate  symmetry  was  a  valid  assumption.  In 
this  section  we  will  attempt  to  address  that  question. 
Several  authors  C18),  C34)  address  the  issue  and  conclude 
that  for  a  symmetrically  stacked  ply  layup,  as  in  this 
research,  quarter  plate  symmetry  can  be  used  because  the 
coupling  stiffness  terms,  ^ ,  are  all  zero.  While  this  is 
true,  the  bending  stiffness  terms,  ^ ,  are  not  all  zero. 

And  while  the  D^ ^  terms  are  not  included  computationally  due 
to  the  numerical  approach  in  the  present  formulation,  they 
serve  as  a  good  baramoter  of  the  degree  to  which  coupling 
between  bending  and  twisting  will  occur.  The  terms  that 
generate  the  problems  with  reguards  to  symmetry 
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In  order  to  study  this,  the 


considerations  are  and 

worst  case,  i.e.  largest,  ratio  of  and  to  was 

identified  and  a  full  and  quarter  plate  model  for  both  a 
linear  and  nonlinear  analysis  were  tested.  Values  of  both 
A  and  D  .  are  shown  in  tables  3.4  and  3.5  for  both 

>-  j  »•  j 

materials  and  all  ply  layups. 


Table  3.4.  Extensional  stiffesses*^ 


Mat. 

Ply 

A.i 

^12 

^22 

^60 

\o 

^2  6 

A 

[0/90]*" 

8 

3.2451 

.0300 

1.6826 

.  0600 

0.  0 

0.0 

A 

C0/±60]^ 

s 

2.6159 

.4639 

1.4441 

.  4938 

0.0 

0.  0 

A 

[ 0/±45 ] ^ 

8 

3.0572 

.6085 

.  7136 

.  6385 

0.0 

0,0 

A 

[0/±45/90]“ 

8 

2.8594 

.4157 

1.2970 

.  4456 

0.0 

0.0 

B 

[ 0/90 ] 

8 

3.3240 

.0901 

1.8229 

.  1920 

0 

0 

0 

0 

B 

[0/±60] 

8 

2.7428 

.4837 

1.6169 

.5856 

0 

0 

0 

0 

B 

t0/±45] 

S 

3.1745 

.6149 

.9227 

.7168 

0 

0 

0.0 

B 

[0/±45/90] 

8 

2.9741 

.  4400 

1.4730 

.  5419 

0.  0 

0.0 

^all  values  are  10® 

'^actual  ply  layup 
‘^actual  ply  layup  [0^ 

^actual  ply  layup  tO, ^X4S^/-45^3 
“actual  ply  layup  t 0^ 2/45^/— 45^/90^] ^ 


Table  3.5.  Bending  stiffnesses** 


Mat. 

Ply 

*^22 

*^66 

*^16 

*^26 

A 

[0/90]"" 

s 

22243 

144 

1410 

288 

0 

0 

A 

[0/±60]^ 

• 

20446 

665 

2166 

809 

231 

682 

A 

[0/±45]“ 

8 

20976 

838 

1289 

982 

527 

527 

A 

[0/±45/'90]‘* 

8 

20949 

813 

1367 

956 

312 

312 

B 

[0/90] 

8 

22361 

432 

2345 

922 

0 

0 

B 

[0/±60] 

8 

20662 

905 

3098 

1394 

234 

643 

B 

[0y±45] 

8 

21180 

1062 

2265 

1551 

507 

507 

B 

[0/±45/90] 

8 

21154 

1039 

2339 

1528 

300 

300 

^actual  ply  layup 
*^act,ual  ply  layup  [0^  ^/60./-60_] 
^actual  ply  layup  [0 ^ ^^45^-45^1 
“actual  ply  layup  [0^_/45  /-45^/90] 

r  j  j  t-  4.  4  a 


•  As  can  be  seen  from  table  3.5,  the  worst  case  was  chosen  as 
material  A  with  a  ^0 ^^y60^y-60^] ^  ply  layup.  In  the  test 
case,  the  plate  was  simply  supported,  with  an  S®30,  and  an 
aspect  ratio  of  one.  The  full  plate  was  modeled  with  an  8x8 
mesh  and  the  1/4  plate  was  modeled  with  a  4x4  mesh.  Both 
linear  and  nonlinear  analysis  were  conducted  and 
nondimensionalized  center  transverse  deflections  were 

#  compared.  The  deflections  were  nondimensionalized  according 
to  Eqn  C3.6>  where  w^  represents  the  deflection  at  the 
center  of  the  plate.  Results  are  shown  in  Table  3.6. 


w  h®E, 


X  100 


C3. 6) 
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Table  3.6.  Center  deflection  results  of  quarter  plate 
symmetry  verification  tests. 


w  Cl inear > 

w  Cnonlinear> 

Full  Plate 

.4573 

.2275 

1/4  Plate 

.4856 

.  2334 

a  Difference 

5.82 

2.  55 

As  can  be  seen  from  Table  3.6,  there  is  i:ood  agreement 
between  the  the  full  and  quarter  plate.  The  author  feels 
that  quarter  plate  symmetry  is  an  accurate  way  to  model  the 
plate  based  on  the  following  factors:  all  comparisions  of 
results  will  be  made  within  this  thesis  and  not  to  other 
published  sources  Cprimarily  because  there  aren't  any  other 
published  cases  for  the  parameters  studied>,  and,  the 
primary  emphasis  in  this  thesis  is  on  trends  in  the  results 
and  not  absolute  numbers. 


IV.  RESULTS/DISCUSSION 

In  this  chapter  the  main  results  of  the  thesis  are 
presented.  The  effects  of  material  properties,  aspect 
ratio,  and  ply  layup  are  discussed  in  separate  subsections. 
The  effects  of  the  type  of  analysis  used,  plate  thickness, 
and  boundary  conditions  lend  themselves  to  being  discussed 
as  they  come  into  play  with  respect  to  the  parameter  being 
analyzed.  Most  of  the  results  in  this  thesis  are  presented 
in  the  form  of  graphs  of  the  nondimensionalized  center 
deflection,  w,  versus  thickness  ratio,  S.  The  center 
deflection  w  is  nondimensionalized  according  to  Eqn  C4.1); 

w  h^E^ 

w  “  — 2 - jJ— ,  X  10  C4.1> 

There  are  three  types  of  forces  that  characterize  the 
response  of  the  plate;  membrane  or  inplane,  bending,  and 
shear.  In  a  linear  analysis  the  response  is  characterized 
by  only  bending  and  shear.  This  is  because  the  membrane 
forces  are  a  function  of  the  membrane  strains  which  are 
neglected  in  a  linear  analysis  due  to  the  transverse 
loading.  In  the  nonlinear  analysis  all  three  types  of 
forces  are  present.  However,  a  thick  plate  is  charaterized 
primarily  by  bending  and  shear  because  there  is  very  little 
displacement  in  the  inplane  directions,  making  the  membrane 
forces  small  compared  to  bending  and  shear.  This  makes  it 
similar  to  the  linear  analysis.  A  thin  plate  is 
characterized  primarily  by  membrane  and  bending  forces. 
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Effects  of  Material  Properties. 


The  graphs  discussed  here  are  representative  of  all  the 
cases  studied.  For  a  complete  set  of  graphs  for  w  versus  S 
considering  all  ply  layups  and  aspect  ratios  that  show  the 
difference  between  the  two  materials  see  Appendix  A.  Only 
one  case  is  discussed  here  for  each  boundary  condition 
because  the  volume  of  data  would  be  overwhelming  and  not 
convey  much  more  information  as  one  representative  graph. 

Simply  supported  plate.  Figure  4.1  is  a  plot  of  w 
versus  S  for  a  simply  supported  plate  with  a 
[0^^/45^/-45^/90^] ^  ply  layup  and  an  aspect  ratio  of  1/2. 
Linear  and  nonlinear  results  for  both  materials  A  and  B  are 
presented.  The  linear  solution,  i.e.  solid  line, 
asymptotically  approaches  the  classical  laminated  plate 
theory  CCLPT)  solution,  which,  as  will  be  shown  later,  is 
constant  for  all  values  of  S,  as  the  plate  gets  thinner. 

This  is  due  to  the  fact  that,  as  the  plate  gets  thinner,  the 
bending  moments  begin  to  dominate  the  response  of  the 
structure  and  the  solution  degenerates  to  a  CLPT  problem. 
Material  A  approaches  the  CLPT  solution  at  approximately 
S=40,  while  material  B  approaches  the  CLPT  solution  at 
approximately  S=25.  As  can  be  observed,  as  the  plate  gets 
thicker,  the  nondimensionalized  center  deflection  increases 
for  both  materials  from  the  CLPT  solution  because  the 
transverse  shear  effects  are  coming  more  and  more  into  play. 
Therefore,  the  main  difference  between  the  linear  solution 
and  the  CLPT  solution  is  that  the  effects  of  the  transverse 
shear  strains  have  been  included.  Note  that  for  any  given 
value  of  S,  material  A  has  a  greater  deflection  than 
material  B.  This  is  because  B  has  a  higher  transverse  shear 
strength,  i.e.  lower  ratio,  than  A.  Consequently, 
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ooooo  A,  LINEAR 
oq_oooA,  NONLINEAR 
B,  LINEAR 

AA_A  A_A  B,  NONLINEAR 


FIGURE  4,1.  w  vs  s  for  a  simply  supported 

[0 ^ _/4S_/-45.]  plate  with  an  aspect  ratio  of  1/2 


when  the  transverse  shear  has  a  major  impact  on  the  response 
of  the  plate,  the  ratio  difference  between  the  two  materials 
becomes  more  important.  Note  also  that  as  the  plate  gets 
thinner,  the  two  linear  solutions  approach  the  same  value. 
Again,  this  is  because  the  bending  forces  become  dominant 
and  the  major  difference  between  the  two  materials,  the 
Ej/Qjg  ratio,  becomes  less  important  as  the  magnitude  of  the 
transverse  shear  forces  decrease  relative  to  the  membrane 
forces.  However,  there  is  always  some  difference  between 
the  solutions  because  there  are  differences  in  the  smaller 
values  of  D  ,  as  can  be  seen  in  table  3.4. 

t  j 

The  nonlinear  solutions,  i.e.  dashed  lines,  also 
asymptotically  approach  a  constant  value,  but,  this  constant 
value  is  not  reached  until  the  plate  is  quite  thin,  i.e.  S 
is  greater  than  60.  This  occurs  when  the  higher  order 
inplane  strains  are  transformed  into  axial  forces.  The 
reason  the  deflections  in  the  nonlinear  solution  are  so  much 
less  than  the  linear  solutions  is  because  the  higher  order 
terms  present  in  the  inplane  deflections  translate  into 
increased  membrane  stiffness  which  plays  a  greater  role  as 
the  plate  gets  thinner  and  the  inplane  stresses  begin  to 
dominate  the  structural  response.  Note  however,  that  for 
values  of  S  between  10  and  20,  there  is  virtually  no 
difference  between  the  linear  and  nonlinear  solution.  This 
is  because  in  this  range,  the  transverse  shear  strains, 
which  are  linear  for  both  geometrically  linear  and  nonlinear 
solutions  Csee  Eqn  C2.10)>,  dominate  the  response  of  the 
structure.  From  the  plot  of  the  curve  in  this  range,  the 
nonlinear  solution  would  appear  to  have  a  greater  deflection 
than  the  linear  solution.  This  is  not  actually  the  case 
however,  and  this  is  a  result  of  the  functions  that  had  to 
be  used  to  plot  the  various  curves.  If  the  actual  data 


points  are  compared,  it  can  be  seen  that  for  small  values  of 
S,  the  deflections  of  the  nonlinear  case  are  equal  to  that 
of  the  linear  case  and  for  larger  values  of  S,  the  nonlinear 
is  less  than  the  linear.  All  simply  supported  plates 
investigated  exhibited  similar  behavior. 

Clamped  Plate.  Figure  4.2  is  a  plot  of  w  versus  S  for 
a  clamped  plate  with  a  f layup  and  an 
aspect  ratio  of  one.  Many  of  the  same  things  can  be  said 
about  fig  4.2  that  were  said  about  fig  4.1,  i.e. 
assymptotically  approaching  the  CLPT  solution,  effect  of 
transverse  shear,  etc.  Here  it  can  be  seen  that  in  the 
linear  solutions,  material  A  is  close  to  the  CLPT  solution 
for  values  of  S  greater  than  SO  and  material  B  for  values  of 
S  greater  than  40.  In  the  nonlinear  case,  the  clamped  plate 
is  similar  to  the  simply  supported  plate,  in  that  the 
solution  does  not  approach  a  constant  value  until  the  plate 
is  quite  thin,  i.e.  for  values  of  S  greater  than  60. 

However,  the  most  striking  thing  about  this  graph  is  the 
relative  closeness  of  the  linear  and  nonlinear  solution  as 
compared  to  the  simply  supported  case  of  fig  4.1.  This 
closeness  is  due  to  the  significant  stiffening  which  the 
clamped  boundary  condition  imparts  to  the  plate. 
Consequently,  the  higher  order  terms  do  not  play  as  large  a 
role  in  the  deflections  and  ensuing  strains  and  stresses. 


Figure  4.3  is  a  plot  of  w  versus  q  for  a  simply 
supported  plate  with  a  C0^^/45^X-4S^/90^3^  ply  layup  and  and 
aspect  ratio  of  one.  The  nondimensionalized  incremental 
load,  q,  is  calculated  by  Eqn  (4.2>  where  q^  represents  the 
load  in  any  given  increment: 


C4. 2) 


For  a  complete  set  of  w  vs  q  graphs,  see  Appendix  B.  The 
results  presented  come  from  the  nonlinear  solution  because 
it  is  Incremental  in  nature.  Solving  Eqn  C4.1)  for  w^  in 
terms  of  w,  one  obtains 


w 

c 


a 


4- 


lOh^E^ 


C4. 3) 


As  can  be  seen  from  Eqn  C4.3>,  the  value  of  w  represents  the 
slope  of  the  actual  load  versus  displacement  curve  because 
E^  and  h  are  constant  and  a  is  constant  for  a  given  value  of 
S.  Based  on  fig  4.4,  which  shows  the  actual 
load-displacement  curve  for  the  same  plate,  it  can  be  seen 
that  as  the  slope  of  the  curve  gets  more  nonlinear  and 
flatter,  the  plate  gets  thinner.  This  occurs  because  the 
higher  order  terms  in  the  assumed  displacement  function  are 
coming  more  into  play  as  the  membrane  forces  begin  to 
dominate  the  solution.  Again,  considering  figure  4.3  for  a 
thick  plate,  i.e.  S=10,  w  is  nearly  constant.  This  implies 
that  w^ ,  <figure  4.4?  the  actual  displacement,  increases 


linearly  with  the  incremental  load,  q^.  This  is  because  for 
thick  plates  the  linear  and  nonlinear  solutions  are 
identical  which  can  be  seen  in  figures  4.1  and  4.2.  Note 
that  as  the  plate  becomes  thinner,  the  difference  in  the 
curves  for  different  material  properties  is  reduced.  This 
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is  because  as  the  plate  gets  thinner  the  membrane  forces 
have  more  effect  than  the  transverse  shear  forces  on  the 
deflections  of  the  plate  and  therefore,  the  major  difference 
between  the  materials,  the  ratio,  has  less  effect. 

Effect  of  Aspect  Ratio. 

Next,  we  look  at  the  effect  of  aspect  ratio  on  the 
response  of  the  plate.  To  serve  as  a  relative  basis  of 
comparison  between  the  various  aspect  ratios,  a  classical 
plate  theory  solution  was  conducted  for  a  simply  supported 
isotropic  plate.  Deflections  were  calculated  using  Navier’s 
solution  by  double  trigonometric  series,  shown  in  Eqn  C4.4>, 
as  presented  in  Szilard  (44): 

03  00 

w(x,y)  =  2  2 

Dn  m  =  t  n  =  1 
odd  odd 

where  m,n  «  1,3,5,. ..,101  and 

Eh® 

D  =  - - -  (4.4a) 

12(l-v^) 

and  for  a  uniform  load. 


2v  2 


mn 


|J!L.  4.  JSlJ. 
^  a^  b^J 


mnx  nrry 

sin  ^  sin 


(4.  4) 


P 

mn 


16q 

re^mn 


(4.4b) 


The  material  properties  chosen  are  shown  in  Eqn  (4.5); 


E  =  20  X  10°  psi,  h  =  . 24  in. 


(4 . 5) 


=  .25, 


q^  =  500  psi 


the  nondimensional ized  deflection  was  calculated  by  Eqn 
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<4.1>  and  is  shown  in  table  4.1. 


Table  4.1.  Nondimensionalized  center  deflection”^ 
Cclassical  solution)  of  a  simply  supported  isotropic 


s 

wCAR*.5) 

w<AR»l> 

w CAR-2) 

10 

.  7851 

.  2968 

.  0491 

20 

.  7851 

.2968 

.0491 

30 

.  7851 

.2968 

.  0491 

40 

.  7851 

.2968 

.  0491 

50 

.  7851 

.  2968 

.  0491 

60 

.  7851 

.2968 

.  0491 

all  w  values  are  times  10' 


As  can  be  seen  in  table  4.1,  as  the  aspect  ratio  increases, 
w  decreases  but,  as  mentioned  before,  it  is  the  same  for  all 
thickness  ratios.  The  value  of  w  for  an  aspect  ratio  of  1/2 
is  roughly  three  times  of  that  for  an  aspect  ratio  of  1.0. 

In  addition,  w  for  an  aspect  ratio  of  two  is  roughly  1/S  of 
that  for  an  aspect  ratio  of  one.  Thus,  one  can  see  the 
significant  differences  aspect  ratio  makes  for  an  isotropic 
plate  of  the  dimensions  given  in  table  3.3. 

The  aspect  ratio  affect  reported  on  next  is  observed 
for  both  materials  and  all  ply  layups.  A  complete  set  of 
plots  for  both  materials  and  all  aspect  ratios  can  be  seen 
in  Appendix  C.  A  representative  graph  for  both  simply 
supported  and  clamped  boundary  conditions  are  discussed 
next. 

Simple  Supports.  Figure  4.5  shows  w  versus  S  for  a 
simply  supported  plate  with  a  layup  made  of 

material  B.  The  curves  for  both  the  linear  and  nonlinear 
solutions  for  all  three  aspect  ratios  are  shown.  Like 
figure  4.1,  the  linear  and  nonlinear  solutions  are  very 
close  for  S  values  less  than  15  to  20.  However,  the  most 
striking  thing  about  this  plot  is  the  proximity  of  the 
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solutions  for  aspect  ratios  of  one  half  and  one  using  linear 
theory  when  compared  to  the  relative  magnitudes  of  w  for  an 
isotropic  plate  as  shown  in  table  4.1.  It  can  be  seen  for 
S=40,  the  value  of  w  for  an  aspect  ratio  of  1/2  is  109/i  of  w 
for  an  aspect  ratio  of  one,  and  w  for  an  aspect  ratio  of  two 
is  38fi  of  w  for  an  aspect  ratio  of  one.  This  similarity  of 
results  between  the  AR=l/2  and  AR»1  cases  occurs  because, 
for  these  two  aspect  ratios,  there  are  a  greater  number  of 
fibers  in  the  direction  of  the  short  dimension  of  the  plate. 
Referring  to  the  discussion  of  ply  layups,  at  least  half  of 
the  fibers  in  all  ply  layups  is  in  the  x-direction,  and 
referring  to  table  3.3,  the  x-direction  is  the  shorter 
dimension  for  an  aspect  one  half  and  is  the  longer  dimension 
for  an  aspect  ratio  of  two.  The  short  dimension  of  the 
plate  has  a  greater  influence  on  the  response  of  the  plate 
and  so,  when  the  short  dimension  of  an  orthotropic  plate 
with  an  aspect  ratio  of  one  half  is  equal  to  the  dimensions 
of  a  square,  i.e,  AR=1 ,  orthotropic  plate  of  the  same  ply 
layup,  one  would  expect  their  response  to  be  very  similar. 
The  nonlinear  solutions  converge  to  values  of  w  that  result 
in  ratios  that  are  somewhat  closer  to  the  isotropic  case 
than  the  linear  solution.  Here,  w  for  AR“l/2  is  16iH  of  w 
for  AR=1 ,  and  w  for  AR=2  is  62%  of  w  for  AR*1 .  Again,  this 
occurs  because  the  higher  order  terms  in  the  inplane 
displacement  function  significantly  stiffen  the  plate  and 
transverse  shear  deflection  effect  is  reduced.  Like  the 
isotropic  case,  as  aspect  ratio  increases,  w  decreases.  All 
simply  supported  plates  exhibit  similar  behavior. 

Clamped  Supports.  Figure  4.6  shows  w  versus  S  for  a 
clamped  plate  with  a  [ 0 ^ 2/45^/-45^/90^ ] ^  ply  layup  made  of 
material  B.  The  curves  for  both  the  linear  and  nonlinear 
solutions  for  all  three  aspect  ratios  is  shown.  Again 
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□  DODD  aR=1  ,  LINEAR 
°°_oqj3AR=1,  NONLINEAR 
ooooo  AR=1/2,  linear 
oooo_o  aR=1/2,  NONL.INEAR 
AR  =  2,  LINEAR 

aa_a^_aar=2,  nonlinear 


FIGURE  4.6.  w  vs  s  for  a  clamped  [0  /45  /-45  /90  1 

plate  of  material  B. 


I 


o 


looking  at  S»40  for  the  linear  results,  the  value  of  w  for 
AR*l/2  is  103k  of  if  for  AR»1 ,  and  w  for  AR"2  is  37k  of  w  for 
AR*1 .  This  compares  favorably  to  the  simply  supported  case. 
The  behavior  here,  in  general,  is  very  similar  to  the  simply 
supported  plate  shown  in  figure  4.5,  but  the  main  difference 
here  is  that  the  linear  and  nonlinear  solutions  are 
comparible  to  a  much  greater  value  of  S.  For  an  aspect 
ratio  of  one  half  the  solutions  are  comparible  up  to  S=40, 
for  AR=1  up  to  S=30,  and  for  AR=2  up  to  S=40.  As  mentioned 
in  the  discussion  of  figure  4.2,  this  occurs  because  the 
boundary  conditions  considerably  stiffen  the  plate  and  the 
higher  order  terms  present  in  the  nonlinear  solution  do  not 
come  into  play  to  the  same  extent  they  do  in  the  simply 
supported  case.  All  the  clamped  plates  studied  exhibited 
similar  behavior. 

Effects  of  Ply  Layup. 

Next,  we  look  at  the  effect  of  ply  layup  on  the 
response  of  the  plate.  Figures  4.7  and  4.8  make  the 
comparison  of  ply  layups  by  plotting  w  versus  S  for  all  ply 
layups  for  both  linear  and  nonlinear  solutions.  The  plate 
shown  is  of  material  A  with  an  aspect  ratio  of  one  and  the 
ply  layups  are  numbered  as  shown  in  table  4.2. 


Table  4.2  Layup  Numbering  System 
for  Figures  4.7  and  4.8. 


Ply 

Layup 

1 

[0.  /90„] 

1  e  8  s 

2 

[0  ^/60./-60.] 

12  6  Os 

3 

[0  „/45./-45.] 

12  6  e  3 

4 

[0  ^/45^/'-45  /90^] 

12  4-  4-  4-3 

2.75  zi 


2.50  d 
3 


W 


2.25 


q 


2.00  q 

I 

1.75  ^ 

1.50  ^ 

3 

1.25  q 

-t 

1.00  ^ 
0.75 
0.50  ^ 
0.25  ^ 
0.00  4^ 


ooooo  PLY:  1 ,  LINEAR 
oo_oooPLY:  1.  NONLINEAR 
a_Q_o_o_p  PLY:  2,  LINEAR 
o«l°«lpPLY:  2.  NONLINEAR 
PLY:  5,  LINEAR 
a^_aaaPlY:  5.  NONLINEAR 
ooooo  PLY:  4,  LINEAR 
oo_oo_oPLY:  4,  NONLINEAR 


Tin  “I"  nTT-T  r  Ti'rf  r"r  "(“ri  r  fr~r  crTT-T' 

0  10  20  30  40  50  60  70 


FIGURE  4.7.  w  vs  s  for  a  clamped  plate 
of  material  A  with  an  aspect  ratio  of  1. 
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Figure  4.7  is  for  a  clamped  plate  and  figure  4.8  is  for 
a  simply  supprted  plate.  As  can  be  seen  in  both  plots, 
there  is  not  a  large  difference  between  the  layups.  This  is 
due  to  the  fact  that  for  each  layup  at  least  S0>£  of  the 
plies  were  in  the  O"  direction,  i.e.  parrellel  to  the 
x-axis,  and  these  plies  dominate  the  response  of  the 
structure.  However,  there  are  some  differences  between  the 
layups  due  to  the  large  differences  between  some  of  the 
smaller  A  and  D  terms  as  a  result  of  the  different  ply 
orientations.  But,  as  can  be  seen  from  tables  3.4  and  3.S, 
these  terms  are  always  at  least  an  order  of  magnitude  less 
than  A^ ^  or  ^ ,  which  dominate  the  plate- s  response,  so  the 
difference  in  deflection  is  small.  The  differences  between 
ply  layups  is  more  pronounced  for  the  nonlinear  solution 
than  the  linear  solution.  Also,  the  differences  are  more 
pronounced  for  a  thick  plate  than  for  a  thin  plate. 


► 
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V.  CONCLUSIONS 


In  this  thesis  the  effects  of  material  properties, 
geometric  nonlinearities,  thickness,  aspect  ratio,  boundary 
conditions,  and  ply  layup  on  the  static  deflections  of  a 
uniformly  loaded  rectangular  composite  plate  were 
investigated.  Parabolic  transverse  shear  strain  was 
included  in  the  analysis.  The  plate  was  analyzed  with  an 
existing  finite  element  code  incorporating  a  4-noded  28  DOF 
rectangular  element.  Quarter  plate  symmetry  was  used  to 
model  the  plate. 

1.  The  linear  solution  asymptotically  approaches  the 
classical  laminated  plate  theory  CCLPT)  solution  as  the 
plate  gets  thinner.  This  is  because  as  the  plate  gets 
thinner  the  inplane  forces  begin  to  dominate  the  response  of 
the  structure  and  the  solution  degenerates  to  a  CLPT 
problem.  As  the  plate  gets  thicker,  the  deflection 
increases  from  the  CLPT  solution  because  the  transverse 
shear  effects  are  coming  more  and  more  into  play. 

Therefore,  the  main  difference  between  the  linear  solution 
and  the  CLPT  solution  is  that  the  effects  of  the  transverse 
shear  strains  have  been  included. 

2.  When  the  transverse  shear  has  a  major  impact  on  the 
response  of  the  plate,  i.e.  the  plate  is  thick,  the 
difference  between  the  two  materials  in  their  transverse 
material  properties  becomes  more  important.  As  the  plate 
gets  thinner,  the  linear  solutions  for  the  different 
materials  approach  the  same  value.  Again,  this  is  due  to 
the  fact  that  the  inplane  forces  become  dominant  and  the 
major  difference  between  the  two  materials,  the 


ratio,  becomes  less  important  as  the  magnitur e  of  the 
transverse  shear  forces  decrease  relative  to  the  inplane 
forces.  However,  there  is  always  some  difference  between 
the  solutions  created  by  the  differences  in  the  other 
material  parameters. 

3.  The  nonlinear  solutions  also  asymptotically  approach 
a  constant  value  but  it  is  much  less  than  the  linear 
solution  and  this  constant  value  is  not  reached  until  the 
plate  is  much  thinner.  The  reason  the  deflections  are  so 
much  less  in  the  nonlinear  solution  than  the  linear 
solutions  is  because  the  higher  order  terms  present  in  the 
inplane  deflections  translate  into  increased  membrane 
stiffness  which  plays  a  greater  role  as  the  plate  gets 
thinner  when  the  inplane  stresses  begin  to  dominate  the 
structural  response.  However,  for  thick  plates  there  is 
virtually  no  difference  between  the  linear  and  nonlinear 
solution.  All  ply  layups  and  aspect  ratios  investigated 
exhibited  similar  behavior. 

4.  For  a  clamped  plate  with  a  linear  solution,  the 
deflection  does  not  approach  a  constant  value  until  the 
plate  is  quite  thin.  However,  the  most  striking  feature 
about  the  clamped  plate  is  the  linear  and  nonlinear  solution 
are  quite  close  compared  to  the  simply  supported  case.  This 
closeness  is  due  to  the  significant  stiffening  which  the 
clamped  boundary  condition  imparts  to  the  plate. 
Consequently,  the  higher  order  terms  do  not  play  as  large  a 
role  in  the  deflections  and  ensuing  strains  and  stresses. 
This  closeness  of  the  linear  and  nonlinear  solutions  was 
exhibited  for  all  ply  layups  and  aspect  ratios  investigated. 

5.  The  slope  of  the  load-displacement  curve  gets  more 
nonlinear  and  flatter,  as  the  plate  gets  thinner.  This 
occurs  because  the  higher  order  terms  in  the  assumed 
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displacement  function  are  coming  more  into  play  as  the 
membrane  forces  begin  to  dominate  the  solution.  However, 
for  a  thick  plate  the  displacement  increases  linearly  with 
load.  This  is  because,  for  thick  plates,  the  linear  and 
nonlinear  solutions  are  nearly  identical. 

6.  For  aspect  ratios  of  one  half  and  one,  the 
deflections  are  quite  similar.  This  is  due  to  the  fact  that 
for  these  two  aspect  ratios,  there  are  a  greater  number  of 
fibers  in  the  direction  of  the  short  dimension  of  the  plate 
for  the  ply  layups  considered  here.  The  short  dimension  of 
the  plate  has  a  greater  influence  on  the  response  of  the 
plate  and  so,  when  the  short  dimension  of  an  orthotropic 
plate  with  an  aspect  ratio  of  one  half  is  equal  to  the 
dimensions  of  a  square  othotropic  plate  of  the  same  ply 
layup,  their  responses  are  very  similar.  The  nonlinear 
solutions  all  converge  to  roughly  the  same  value  for  all 
aspect  ratios.  This  is  because  the  higher  order  terms  in 
the  inplane  displacement  function  significantly  stiffen  the 
plate  and  reduce  transverse  deflection.  For  a  clamped 
plate,  the  linear  and  nonlinear  solutions  are  comparible 
until  the  plate  is  much  thinner  than  for  a  simply  supported 
plate.  This  occurs  because  the  boundary  conditions 
considerably  stiffen  the  plate  and  the  higher  order  terms 
present  in  the  nonlinear  solution  do  not  come  into  play  to 
the  same  extent  they  do  in  the  simply  supported  case.  All 
clamped  plates  exhibit  similar  behavior. 

7.  In  this  research,  ply  layup  does  not  have  a  large 
effect  on  the  deflection  of  the  plate  since  at  least  505<  of 
the  plies  are  in  the  0°  direction,  i.e.  parrellel  to  the 
x-axis,  and  these  plies  dominate  the  response  of  the  plate. 
However,  there  are  small  differences  in  deflection  between 
the  layups.  These  differences  are  more  pronounced  for  the 


nonlinear  solution  than  the  linear  solution  and  are  more 
pronounced  for  a  thick  plate  than  for  a  thin  plate.  Also, 
it  was  observed  that  the  maximum  transverse  shear  stress 
usually  occurred  between  the  plies  with  the  maximum  angle 
difference  between  them. 
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Appendix  A 


Material  Properties  Comparisons 


This  appendix  presents  a  complete  set  of  curves  of 
nondimesionalized  displacement,  W,  versus  thickness  ratio, 

S,  for  all  ply  layups  and  aspect  ratios.  Linear  and 
nonlinear  solutions  for  both  materials  are  presented  on  each 
graph.  The  displacement  is  nondimensional ized  according  to 
Eqn  CA.  1): 


where 


w  h=^E, 

w  =  — S - -i—  X  10  CA.l) 

=  center  displacement 
h  =  plate  thickness 

Ej  =  Young’s  modulus  in  fiber  direction 
=a  total  load  on  plate 
a  =  characteristic  inplane  dimension 
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Figure  A. 3.  w  vs  s  for  a  simply  supported 
[0  /4S  /-45  ]  plate  with  an  aspect  ratio  of  1 
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Figure  A. 7.  w  vs  s  for  a  clamped 

[0  /45  /-45  ]  plate  with  an  aspect  ratio  of  1 
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Figure  A. 11.  w  vs  s  for  a  simply  supported 
CO  /43  /-43  3  plate  with  an  aspect  ratio  of  1/2 
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Figure  A. 18.  w  vs  s  for  a  simply  supported 
CO^  /60  /-60  1  plate  with  an  aspect  ratio  of  2 
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Figure  A. 22.  w  vs  s  for  a  clamped 

CO  /60  /-60  ]  plate  with  an  aspect  ratio  of  2 

12  6  6  a 


ooooo  A,  LINEAR 
oo_oo_oA,  NONLINEAR 
B,  LINEAR 
B,  NONLINEAR 


Figure  A. 24.  w  vs  s  for  a  clamped 

CO  /45  /-43,/90,]  plate  with  an  aspect  ratio  of  2 


Appendix  B 


Nondimensionalized  Load-Displacefflent  Curves 


This  appendix  presents  a  complete  set  of  curves  of 
nondimesionalized  displacement,  w,  versus  nondimensionalized 
load,  q,  for  all  ply  layups  and  aspect  ratios.  Nonlinear 
solutions  for  both  materials  are  presented  on  each  graph  for 
thickness  ratios  of  10,  30,  and  60.  The  displacement  and 
load  are  nondimensionalized  according  to  Eqn  CB.l): 


where 


w  h  E,  q  _ 

w  »  - 2 - -i_  X  10  q  »  X  10®  CB.l) 

q<,a  1 

=“  center  displacement 
h  *  plate  thickness 

E^  =  Young's  modulus  in  fiber  direction 
q^  *  total  load  on  plate 
a  =  characteristic  inplane  dimension 
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Figure  B.2.  w  vs  q  for  a  simply  supported 
[0  /60  /-60  ]  plate  with  an  aspect  ratio  of  1 
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Figure  B. 6.  w  vs  q  for  a  clamped 
[0  /60_/^-60^1  plate  with  an  aspect  ratio  of  1 
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Figure  B.IO.  w  vs  q  for  a  simply  supported 
[0  /60  /-60  1  plate  with  an  aspect  ratio  of  1/2 
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Figure  B.ll.  w  vs  q  for  a  simply  supported 
[0. „/45./-45.] _  plate  with  an  aspect  ratio  of  1/2 
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Figure  B.12.  w  vs  q  for  a  simply  supported 

[0  /45  /-43  /90,1  plate  with  an  aspect  ratio  of  l/'2 
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Figure  B. 13.  w  vs  q  Tor  a  clamped 

CO  /90  ]  plate  with  an  aspect  ratio  of  1/2 


Figure  B.14.  w  vs  q  for  a  clamped 

[0  /60  /-60  J  plate  with  an  aspect  ratio  of  1/2 
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Figure  B.13.  w  vs  q  Tor  a  clamped 

[0  /45  /-45  ]  plate  with  an  aspect  ratio  of  1/2 
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Figure  B.16.  w  vs  q  for  a  clamped  ^ 

[0  /45  /-45  /90  ]  plate  with  an  aspect  ratio  of 
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Figure  B.17.  w  vs  q  for  a  simply  supported 
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Figure  B.18.  w  vs  q  for  a  simply  supported 
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Figure  B.19.  w  vs  q  for  a  simply  supported 
CO  ^45  /-45  ]  plate  with  an  aspect  ratio  of  2 
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Figure  B.23,  w  vs  q  for  a  clamped 

[0  /43  /-45  ]  plate  with  an  aspect  ratio  or  2 
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Appendix  C 


Aspect  Ratio  Comparisons 


This  appendix  presents  a  complete  set  of  curves  of 
nondimesionalized  displacement,  W,  versus  thickness  ratio, 
S,  for  all  ply  layups.  Linear  and  nonlinear  solutions  for 
each  aspect  ratio  are  presented  on  each  graph.  The 
displacement  is  nondimensionalized  according  to  Eqn  CC.l): 


where 


w  h^E, 

w  »  — 2 - L_  X  10  CG.l) 

w^  “  center  displacement 
h  =  plate  thickness 

Ej  *  Young's  modulus  in  fiber  direction 
=  total  load  on  plate 
a  *  characteristic  inplane  dimension 
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Figure  C.3.  w  vs  s  for  a  simply  supported 
[0  /• — 43^1  plate  made  of  material  A 
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Figure  G.4.  w  vs  s  Tor  a  simply  supported 
CO, ]  plate  made  of  material  A 
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Figure  G.8.  w  vs  s  for  a  clamped 
[0  ]  plate  made  of  material  A 
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